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ABSTRACT 

I present exact expressions for the interior gravitational potential \^ of a 
system of concentric const ant- density (Maclaurin) spheroids. I demonstrate 
an iteration procedure to find a self-consistent solution for the shapes of the 
interfaces between spheroids, and for the interior gravitational potential. The 
external free-space potential, expressed as a multipole expansion, emerges as 
part of the self-consistent solution. The procedure is both simpler and more 
precise than perturbation methods. One can choose the distribution and mass 
densities of the concentric spheroids so as to reproduce a prescribed barotrope 
to a specified accuracy. I demonstrate the method's efficacy by comparing its 
results with several published test cases. 

Subject headings: Planets and satellites: interiors 



1. Introduction 

In its general form, the problem of the theory of figures is to find the external gravita- 
tional potential of a liquid planet in hydrostatic equilibrium, rotating at a uniform rate u, 
and obeying a specified barotropic relationship for the dependence of pressure P on mass 
density p. 

The expected precision (~ one part in 10^) of the Juno Jupiter orbiter spacecraft's 
measurements of Jupi ter's gravity field will require a gra vitational-modeling theory of un- 



precedented accuracy (jKaspi et al.ll2010l ). iHubbardI (l2012l ) (Paper I) is an initial step toward 



such a theory. Paper I presents a new approach to the calculation of the multipole expansion 
of the external gravitational potential of a rotating planet in hydrostatic equilibrium. 

As is well known, the problem of the theory of figures can be solved in closed or 
partially-closed form for a small number of special barotropes but arbitrary barotropes 
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generally require numerical methods. Analytic methods balloon in complexity even for 
the apparentl y simple case of two constant-density layers, the so-called two- layer Maclau- 



rin sp heroid (ISchubert. Anderson. Zhang. Kong, fc Helledl 1201 ll : iKong. Zhang, fc Schubert 



20121 ). In principle, such analytic complexity could be bypassed by seeking a purely nu- 
merical solution to the general equation of hydrostatic equilibrium. However, numerical 
solutions are vulnerable to numerical noise, produced for example by cancellation of nearly 
equal terms. In the traditional approach, cancellation of terms is mitigated by solving a hi- 
erarchy of integrodifferential equations gene rated from a perturbation e xpansion of the mass 



distribution in powers of the rotation rate ( iZharkov fc Trubitsyrull978l ) 



Paper I shows how the particular problem of the hydrostatic equilibrium of a rotat- 
ing const ant- density planet can be numerically solved to high precision by using gaussian 
quadrature to obtain the mass multipole moments. In this method, the moments are calcu- 
lated by performing one- dimensional integrals over the surface mass distribution. Although 
gaussian quadrature is a numerical approximation to analytic integration, the results are 
exact (to within the floating point precision of the computer), as long as the integrand can 
be expressed as a polynomial of degree less than the degree of the gaussian quadrature; for a 
Maclaurin spheroid, using this approach with 48 quadrature points yields numerical results 
with a precision of at least ~ 10^^^. The mass multipole moments are then self-consistently 
iterated on the shape of the surface (Paper I). The method of Paper I largely bypasses 
the analytic complexity of perturbation methods and avoids cancellation problems, but as 
presented is only valid for a constant-density object, or for a const ant- density object with 
special boundary conditions. 

The present paper shows how the method of Paper I is straightforwardly generalized 
to solve the problem of multiple-layered constant-density spheroids. The resulting method, 
which I call the concentric Maclaurin spheroid (CMS) method, retains all of the advan- 
tages of the approach of Paper I, with the additional flexibility that the concentric Maclau- 
rin spheroids can be arranged in sufficient numbers to closely approximate any prescribed 
barotrope. As we will see, an actual density discontinuity such as a discrete core or first- 
order phase transition is trivially incorporated in the CMS method, as opposed to the usual 
theories of figures. 

In the following Section 2, I present the analytic development of the CMS method. In 
Section 3, I apply the CMS method to several published test cases and I show how the 
method can incorporate a prescribed barotrope. In the conclusion (Section 4), I discuss 
how the CMS method can be applied to analysis of Juno gravity data expected to arrive 
beginning in 2016. 
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2. Theory for N Layers of Maclaurin Spheroids 



2.1. Exact calculation of gravitational potential 



Consider a configuration of N concentric Maclaurin spheroids (Fig. 1). Label the 
spheroids with index i = 0, 1,...,A^ — 1, with i = corresponding to the outermost spheroid 
and i — N — 1 corresponding to the innermost. 

Because the gravitational potential V is linear in the mass density p, we may use the 
principle of superposition, such that the total potential at any point in space is the sum of 
the partial potentials of N concentric constant-density spheroids. Figure 2 illustrates this 
concept for a three-layer model. 

Let the equatorial radius of the outermost spheroid be Oq, and let the equatorial radii 
of the concentric spheroids be ao > ai > 02 > . . . > ajv-i. The total external gravitational 
potential at some point "A" on the outermost level surface is 



k=0 / 

where r is the radius from the center of the planet, /i is the cosine of the angle from the 
rotation axis, the P2fc(A*) a-re the usual Legendre polynomials. 



etc., where the relation = rj(/x) is the surface equipotential of the i-th layer. The zero- 
degree values are given by 




(1) 






(3) 




(4) 




(5) 



and so we have for the total mass M 
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Fig. 2. — Method of superposition of Maclaurin spheroids, for the case N — 3. The point 
"A" is a typical point on the outermost surface of the planet. 
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N-l 



M = Y, A,o. (6) 



i=0 

We now introduce the usual dimensionless forms of the multipole moments, 

MafJi,2k = -A,2fe. (7) 
and the dimensionless radii of level surfaces 

= ri{l^)/ao- (8) 
The total external gravitational potential at point "A" can thus be rewritten 

Kxt,A = ^ 1 - E E -^^.^^ M'^'p^M ' (9) 



=0 k=l 



where 



for i > and 



Spi ^ Pi- Pi-i (11) 



5po = po- (12) 



Next, wc must compute the total gravitational potential on an interior interface (level 
surface) at an arbitrary point "B", as shown in Fig. 3. 

First, we consider a problem in which there is only a mass distribution with a constant 
density 6pi interior to point "B" located at coordinates {r,p). We calculate the external 
potential due to this mass distribution, finding 

„ oo 

V^,e.t,B = -Yl D^,2kr-'''P2k{p). (13) 
k=0 

Next wc calculate the external potential at the surface of a spherical mass distribution with 
radius r and constant density (shown as a dashed circle in Fig. 3): 

^i,ext,B ~ ^ 3 ^ ^ 
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Fig. 3. — Schematic diagram illustrating the computation of three contributions to the 
gravitational potential at point "B" on an interior interface. 
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Finally, we calculate the internal potential at point "B" due to the mass distribution with 
constant density external to the dashed circle in Fig. 3: 



V'. , o = 



/ df,'P2,{pi') dr'r'-''^'. (15) 



Adding all contributions to the potential at point "B" due to the mass density in the 
i-th layer and in the i — 1-th layer, one has 

^ oo oo 

^^.^ = 7 E D^,2kr-''P2M + Gj2DU,2kr''P2M + GD'U^y, (16) 

fc=0 fe=0 



where [c/Eq. (4)] 



for A; > 1 we have 



ioY k — 1 we have 



and for /c = 



A,2. = ^ c^/. PM r,(/.)^^+^ (17) 







^:-i,2. = rf/^ i^2fe(/.) r._l(/.)^-^^ (18) 



2 



= 47r(5p,_i / ci//P2(//)ln[ri_i(//)], (19) 
D^^o = 27r5pi_i /" rf/iri_l(/i)^ (20) 

^0 



d:u,„=-^. (21) 

Next, analogous to Eq. (7), we introduce dimensionless forms of the -D^_i^2fe- 

M%'-^^Jl,, = -Dl^,„ (22) 

and 

By analogy with Eq. (10), we may write the dimensionless forms of Eqs. (18-21): for 
A: > 1 

ior k — 1 

J, ^ I ' ^Pi lo dfJ' f2(/^) HM] \ 
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and for A; = 



J' 



7" 







(26) 

(27) 



Eq. (16) then takes the form 
GMl 



.k=0 



k=0 



(28) 



The total potential at a point B located at coordinates (^, /j.) on the j'-th interface is 
obtained by summing over all layers: 



VbU) = 



GMl 



do ^ 



N-l 00 

$]$^W-''^2.(/^) 

i=j k=0 



i-i 00 

+ E E ^;-i.2.e^'=^'^2.(/.) + E ^^-1/ 

i=0 fe=0 i=0 



(29) 



2.2. Parameters and scaling 

Assume that the planet rotates as a solid body at an angular rate ou. Therefore in the 
corotating frame there appears a rotational potential 

g = irV[l-P2(/^)], (30) 
and the total potential U appearing in the equation of hydrostatic equilibrium 

VP = pVU (31) 

is given by 

U^V + Q. (32) 

For a nonrotating planet, all multipole moments for A; > vanish and the potential V within 
the planet depends only on r. The presence of the nonspherical term Q in U breaks the 
spherical symmetry and excites all of the k > terms. We represent the magnitude of Q by 
the dimensionless parameter 
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2^3 

(33) 



'GM 



The number and location of the concentric Maclaurin spheroids can be chosen arbitrarily. 
Let the equatorial radius of the i-th spheroid be a^. Let 

A, = ^,i = 0,l,...,7V-l. (34) 

The Aj can be spaced equally between and 1, or could be made denser in certain regions 
(for example, one could space them at two or three per density scale height). 

Define the mean density of the planet: 

_ 3M 1 

4^«o/o'rf/xCo(/x)3 

For numerical convenience one may use the dimensionless density increment for the i-th 
spheroid: 

5i = Spi/p (36) 

As can be seen by examining Eqs. (10, 24-27), the dimensionless multipole moments 
can be calculated using either the 6pi or the 6i. However, although the moments are di- 
mensionless, further scaling is necessary in order to achieve satisfactory numerical accuracy. 
Consider, for example, a model with N — 128, having spheroids with equally-spaced equato- 
rial radii. It then becomes necessary to consider spheroids with Aj ~ 1/100, so for example 
J{qo2o has an integrand ~ 10"^^^*^"^^^. The resulting huge number is then multiplied by 
~ 10"2^(+^^) in the corresponding term in Eq. (28). To avoid pointlessly multiplying and 
then dividing by large factors, we rescale to new variables and parameters: 

Ciipi) = (37) 



Ji,2k = Ji,2kl^^, (38) 
Ji,2k = Ji,2k^T^^^ (39) 

Then 
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for A; > 1 



for A; = 1 



and for A; = 



[2-2k)[ Ejro^^AVo^/^GW '' ^ ^ 



T^fJo ^1 dfi Cj W 



3 / 5i\ljld^iP2k{^i)Q{^^? ' 

'2 y>:U^^>^'t.d^,Q{^,f 



We introduce dimensionless planetary units of pressure (-Ppu), density (ppu) , and total 
potential (f/pu), such that 

P ^ ^P.. (44) 



«5 



P = —Ppu, (45) 

[/ = ^[/p,. (46) 
ao 

Evaluating the total potential at the surface of the outermost Maclaurin spheroid at the 
equator (/x = 0), we have 



N-l oo 

t^o.pu = 1 + 2? - E E •^^,2.Af P2.(0). (47) 



1=0 k=l 

At the surface of each subsequent Maclaurin spheroid we have 



^ /AT-l oo 

^ \ i=j k=0 

i-1 oo ^ j-1 \ 

+ E E 42.(A./A.)^'^-^^P2.(0) + E J^^>^ + (48) 

1=0 k=0 1=0 / 

and at the center of the planet 

N-l N~l 

^center,pu = ~ '^i,2k^i ~ ~ •^i,2k- (49) 

i=0 1=0 
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The shape Co(a*) of the surface of the planet is an equipotential given by the solution to 

/ iV-l oo \ 
- 1 - E E J^,2k>^^VP2M + ^gCo'[l - P2{fl)] = t/o.pu- (50) 



i=0 k=l 



Correspondingly, the shape Cjil^) of the surface of the j-th spheroid is an equipotential 
given by the solution to 



^N—l oo j — l oo 



^ E E JiAK/Xjrci''P2M + E E JU(^^/^^r^'c?^'pM 

\ i=j k=0 1=0 k=0 

j-1 \ /N-1 oo 

+ E J^,0^<! + - + E E J^M^^/^JrP2k{0) 



i=0 / \ i=j k=0 

j-1 oo j-1 

+ E E J^,2k{^J/^^r'^'P2k{o) + j::oA| ) - ^^aj = o. (si) 

i=0 k=0 i=0 



2.3. Gaussian quadrature 

All of the foregoing expressions for the potential of N concentric Maclaurin spheroids 
are exact. For practical applications, we are interested in finding the potential as a multipole 
expansion to finite (say, thirtieth) degree, corresponding to an upper limit at, say, k^ax — 
15. For this purpose one may numerically evaluate the angular integrals for the multipole 
moments using L > 2kmax gaussian quadrature points. For the examples presented in this 
paper, we use L = 48 gaussian quadrature points fia, oc — 1,2,...L with corresponding 
weights Wa, a — 1,2, ... L over the interval < /i < 1. 

Using initial guesses for the moments Ji,2fe, <^i,2fe) ^'/fij solves Eqns. (50) and 
(51) for C«(A*a) for i = 0, 1, . . . , iV — 1 and a = 1, 2, . . . , L. The solutions for these values are 
then used to evaluate the gravitational moments by gaussian quadrature. 



^i^f Zla=l ^« ^2fc(/^a) (i{^J'a 



,2k+3 



^2k + 3j[ Eto^<^.A3ELi^«0(/^a)^^ ^ ^ 



etc. 



One then iterates between calculation of the level surface shapes via Eqns. (50) and 
(51) and the gravitational moments via Eqns. (40-43) until the difference between successive 
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iterations falls below a specified tolerance. For the purposes of achieving Juno-level precision, 
about 30 such iterations (over all N spheroids) usually suffices. 



2.4. Calculation of the barotrope 

First, we calculate the density in each uniform layer; for the j-th layer 



Efc=0 ^fc-^fcEa=lCfe(/i 



Next, we calculate the total potential t/pu on the outer surface, on each of the interfaces, 
and at the center, using Eqs. (47-49). Since the density is constant between interfaces, Eq. 
(31) is trivially integrated to obtain the pressure at the bottom of the j-th layer: 

Pj,pu Pj~l,pu ~l~ Pj — l,puiUj,pu f^j — l,pu)- (^'^) 

Figure 4 shows an example of the resulting stair-step barotrope obtained for a rotating 
Jupiter model with = 32 and a linear variation of density with mean radius (the linear- 
density model is discussed further below). 



3. Comparison of CMS results with test cases 
3.1. Linear density profile 



Results for a linear density naodel o f Jupiter using a fifth-order theory are tabulated in 
Table 3.1 of IZharkov fc Trubitsynl (119781 ). They adopt a mass-density profile which is linear 
in the mean radius of a level surface rather than in its equatorial radius. The mean radius 
Sj of the j-th level surface relative to the planetary mean radius is given by 



11 
So 



1/3 



(55) 



If we arrange a constant increment 6j in Xj (with constant AA), we can make the resulting 
density linear in s/sq by modifying the density increment for each spheroid to 



5., 



(56) 
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Fig. 4. — Inferred variation of pressure vs. density (both in c.g.s. units) in a CMS model of 
Jupiter with = 32, for an assumed hnear variation of mass density with mean radius. 
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The intervals As must be computed iteratively. Furthermore, IZharkov &: Trubitsynl (119781 ) 
expand their fifth-order theory in powers of the small parameter 



m 



GM' 



(57) 



with a fixed value of m that differs slightly from the value obtained from the value obtained 
for a more realistic Jupiter model. Thus, the CM S calculations must b e also iterated to 
obtain a value for m that matches the one given by lZharkov &: Trubitsynl (1l978l ). 

Table 1 presents a comparison of results for the linear-density model. Agreement is 
excellent for = 128. The inferred pressure- density relation for = 32 was depicted in 
Fig. 4. 



3.2. Two-layer Maclaurin spheroids 

The relative simplicity and elegance of Maclaurin's theory for the single spheroid disap- 
pears for A^ = 2 . Nevertheless , one finds considerable literature for the case N = 2, dating 



back at least to iDarwinI ( 119031 ) 



First, it is useful to test the CMS theory by calculating the equipotential shape of an 
interior interface in a Maclaurin spheroid of uniform density. For this test, I set A^ = 2, 
Ai = 0.5, 6o = 1, 5i = 0. I set q equal to the Jovian value adopted in Paper I. The 
converged CMS model agrees exactly with results presented in Paper I, as it should. Figure 
5 shows the deviations of the outer and intermediate surfaces from ellipsoids of revolution, 
with 6( = C(/i) — 1/ a/(1 + P^i^), where £ is related to m by Maclaurin's result, m = [(3 + 
£^)arctan£ — 3£]. Evidently the shape of the intermediate surface is to high precision an 
ellipsoid of revolution, homologous to the outer surface, as it is in Maclaurin's analytic theory. 



Next, I compare my N = 2 CMS results with the results of lSchubert. Anderson. Zhang. Kong, fc Hellec 



(120111 ). The Schubert et al. models are characterized by three parameters. 



Jo dfi^oifiy 

the core-envelope density ratio pi/po, and a dimensionless rotation parameter 



all in m y notation. In the present paper I compare three models adopted by lSchubert. Anderson. Zhang. Ko 



( I2OIII ): "Mars", "Neptune", and "Uranus2". The quantities that they compute for these 
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Fig. 5. — Departure of outer surface (dashed) and intermediate surface (solid) from an 
ellipsoid of revolution, for a classical Maclaurin spheroid. 
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models are J2, and Eq and Ei, respectively the eccentricities of the outer surface and inter- 
mediate surface, defined by 

E = v/r^"(r^, (60) 

where the oblateness e = 1 — ({fi = 1). The CMS calculations require iteration to match 
the values of and €2- Results are presented in Tables 2-4. While the values for "Mars" 
generally agree, there are unexplained discrepancies for "Neptune" and "Uranus2" . The "3rd 
order" results of Schubert et al. agree with CMS results to high precision, but their "exact" 
results differ by larger-than-expected amounts. 



3.3. Polytrope of index one 

The polytrope of index one is defined by the barotrope 

P = Kp\ (61) 

where the polytropic constant K is chosen in the present application to yield a model planet 
matched to Jupiter's mass and equatoria l radius. Rotating planet models obeying this 



barotrope have been extensively studied (IZharkov fc Trubitsynl Il978l : iHubbardI 119751 ) . so 



it provides a rigorous test of the CMS method. 

Moreover, the study presented in this section provides a useful illustration of how, for 
a chosen barotrope, one may choose CMS arrays of \j and 5j to yield a close match to the 
barotrope. 

For a nonrotating n = 1 polytrope, the density distribution is given by 

sinvrA , , 

P = Pc^-, (62) 

where pc is the central density. To obtain a first approximation to the 5 distribution over 
the spheroids, we differentiate: 

d{p/Pc) costtA sinvrA 

and we use this relation to obtain starting values of the 6j. 

After obtaining a converged hydrostatic-equilibrium model for spheroids with the 
above array of 6j, one calculates the arrays f/j,pu and -Pj,pu- Next one calculates an array of 
desired densities Pj,pu, desired according to 

Pj,pu, desired = P ( ^ + ^i) ) ' (64) 



- 18 - 



where p{P) is the inverse of the adopted barotrope -P(p). Differencing the desired densities 
between layers then gives a new array of 6j. In general, it is necessary to scale the densities 
so as to obtain the correct total mass of the CMS model. This can be effected by rewriting 
the barotrope as 

P = P{Cp), (65) 

where C is a dimensionless factor. For the polytrope of index one, when one adopts a value 
of C greater or less than one, this is equivalent to redefining the value of K. 

After obtaining a new converged CMS model, the process of adjusting the densities to 
obtain a new array of Pj,pu, desired , etc., continues until all changes in gravitational moments 
and in the value of C are reduced to within a specified tolerance. Because of the additional 
step of fitting the barotrope, more iterations are required for convergence. Figure 6 shows 
the fitted and target n = 1 barotrope of a converged 512-layer CMS model of Jupiter. 

The comparis on models for the n = 1 polytrope are (1) analytic expansions of J2, J4, 
and Jq to order (IZharkov fc Trubitsynlll978t lHubbardlll975l ). and (2) a self-consistent-field 
calculation of the rotating polytrope based on the analytic result that the interior density 
can be expaiided a s a series of products of spherical harmonics and spherical Bessel functions 
jn ( lHubbardlll975l ). Table 6 presents results for = 256 and = 512 CMS models along 
with the comparison models. The results for J2 for the A^ = 512 CMS model were still 
changing in the fifth figure after the decimal point after fifteen iterations on the barotrope 
fit. 



3.4. Convergence considerations 



Kong. Zhang, fc SchubertI (120131 ) have criticized the Maclaurin spheroid approach em- 



ployed in this paper, stating that the method of Paper I is incomplete. Further clarification 
is called for since the method of Paper I is central to the CMS method. 



Consider a Maclaurin spheroid of eccentricity 
usual. 



For its external potential write, as 



Vext(?^,/i) 



GM 



k=l 



(66) 



Where does this infinite-series expansion diverge for the Maclaurin spheroid? Evaluate it at 
the spheroid's p ole, where /i = 1 and r^/a^ = b'^/a^ = 1/(1 + i"^). Substitute Eq. (10) of 
Hubbard! J2OI2I 1. We get 
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Fig. 6. — Target (solid line) n = 1 barotrope for a Jupiter model, and fitted = 512 CMS 
barotrope (stairstep). 



-20- 



viM) 



^"^(2fc + l)(2fc + 3) 



6 

The ratio of the k + 1-th to the fc-th term is 

tk+i ^ ^2 (2fc + l)(2A; + 3) 
tfc (2A; + 3)(2A; + 5)' 



(67) 



(68) 



Therefore the series converges if £^ < 1 or the oblateness e < 1 — b/a = 1 — 1/ \/2 



0.29289, in agreement with the estimate of iKong. Zhang, fc SchubertI ( 120131 ). The corre- 



sponding m = merit = 0.212389 and q = gcrit = 0.424778. These values are far larger than 
the parameters of any known planet. Note, by the way, that the point of bifurcation for the 
Maclaurin-Jacobi ellipsoid sequence is at a somewhat larger ^bifurc = 1-39, corresponding to 
"^bifurc = 0.280 and gbifurc = 0.669. 

Paper I shows that for a Maclaurin spheroid with Jupiter's q = 0.089, my method 
gives for the shape of the spheroid's surface a numerical result that differs no more than a 
few parts in 10^'^ from the exact Maclaurin shape. Note that Figure 5 of this paper shows 
similarly-small departures at the outer surface and on an intermediate surface. Thus, for 
this value of q, any neglected terms in Equation (66) will not exceed ~ 10^^^ of the included 
terms. 

Repeating the calculation for a Maclaurin spheroid with a Saturn-like q = 0.155, I 
obtain results shown in Figure 7 for the relative difference of the spheroid's surface radius 
from the exact Maclaurin ellipsoid shape. The departures are, in absolute terms, no more 
than a few centimeters, and would have no significance whatsoever for practical models of 
Saturn's gravity field. Moreover, the real Saturn is much less oblate than the Maclaurin 
model, so the departures of a CMS model from an "exact" model will be smaller still. 

As in Paper I, the general CMS method relies upon the requirement that the exter- 
nal multipole expansion (66) of a given spheroid's potential converges at all points on the 
spheroid's surface. First, on a sphere of radius r = ag, the expansion converges. To see this, 
note that ~ (— l)'^"'"^g'^ for a uniformly-rotating body in hydrostatic equilibrium. Thus, 
on the sphere = 1, the ratio of the k + 1-th term to the k-th term is ~ —q, so for g < 1 
the series converges. 

Next we examine the series convergence at the pole, /i = 1, where = 1 ~ e, where 
e ~ g is the oblateness. At /i = 1, the ratio of the k + 1-th term to the k-th term is 
~ —q/{l — e)^. Thus, as long as g < C"(l — e)^ (where C is a constant of order unity 
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whose precise value depends on the barotrope) the series converges. As discussed, the q 
values for Jupiter and Saturn are such that the convergence criterion is well satisfied (as 
numerically demonstrated fo r the test cases). See also a relevant discussion in Section 38 of 



Zharkov fc Trubitsvnl (119781 ). 



4. Practical application to analysis of gravity data 

The CMS analysis technique presented here can be vectorized for efficiency, although 
no significant effort has been made to do so at this point. For practical computations it will 
probably be necessary to further increase and to increase the number of iterations on the 
barotrope fit, in order to match the theoretical results to the expected precision of spacecraft 
measurements. 

Further iteration loops will be required if a subset of the calculated are to be fitted to 
observed values. Adjustable parameters might include: (1) the mass and density of a discrete 
core at the planet's center, (2) chemical and density discontinuities at various layers, and 
(3) modifications to the assumed barotrope (crudely illustrated in this paper with the scale 
factor C). 

As is obvious, there exists an infinity of possible arrangements of spheroids which can be 
fitted to a finite set of gravity data. Thus, a unique inversion cannot be achieved. However, 
application of specific physically-based barotropes and cosmochemical considerations can 
lead to the most realistic interior models. 



5. Conclusion 



One can further generalize the CMS method in two directions. First, in addition to the 
rotational potential Q one may introduce a tidal potential from a satellite. The resulting tidal 
perturbing potential Qtid will be a function of two angular variables, /x and 0, where is the 
angle from the sub-satellite longitude. Since Qtid will excite both zonal and tesseral gravity 
harmonics, evaluation of the response on all CMS surfaces will require two-dimensional 
rather than one-dimensional integrals. However, there appears to be no practical barrier to 
evaluating such integrals using two-dimensional gaussian quadrature (to be sure, at the cost 
of more computing time). 

Second, one c an investigate the effects of diffe renti al rotation on cy linders (for related 
investigations, see iKong. Zhang, fc SchubertI (120121 ) and iHubbardI ( 1l982l )). 
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Fig. 7. — Calculation of the difference 6^ = ^cms ~ Cexact for a Maclaurin spheroid with 
Saturn's mass, mean density, and rotation rate. 
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Table 1. Linear Density Model 



Quantity 


/ i 5 —order tiieory 


LMb theory (A* = 128) 






0.088822426 


m 


0.0830 


0.082999915 


J2 X 102 


1.4798 


1.4798138 


-J4 X 10^ 


5.929 


5.9269129 


Je X 10^ 


3.497 


3.4935680 


-Jg X 10^ 


2.52 


2.5493209 


Jio X 10^ 


2.4 


2.1308951 


-J12 X 10^ 




1.9564143 


Jl4 X 10^ 




1.9237724 



'Zharkov fc Trubitsvn 



Table 2. "Mars" Note. "3rd order" and "exact" values from Schubert et al. (2013) 



Quantity value "3rd order" "exact" CMS (A^ = 2) 



0.125 

Po/Pi 0.486 

€2 0.00347 

q 0.0046205430 

J2 X 10^ • • • 1823.18 1823.1 1823.1832 

El ••• 0.0888747 0.088859 0.088874693 

Eo ■■■ 0.100295 0.10030 0.10029471 
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Table 3. "Neptune" Note. "3rd order" and "exact" values from Schubert et al. (2013) 



Quantity 


value 


"3rd order" 


"exact" 


CMS (N = 2) 


Qv 


0.091125 








Po/pi 


0.157334 










0.0254179 








q 


0.026207112 








J2 X 10^ 




6188.92 


6241.0 


6188.9267 


El 




0.143515 


0.15147 


0.14351534 






0.209658 


0.21019 


0.20965898 



Table 4. "Uranus2" Note. "3rd order" and "exact" values from Schubert et al. (2013) 



Quantity 


value 


"3rd order" 


"exact" 


CMS {N = 2) 


Qv 


0.0563272 








Pol pi 


0.0791231 










0.0318902 








<1 


0.029581022 








J2 X 10^ 




5680.32 


5801.4 


5680.3242 


El 




0.115655 


0.14160 


0.11565564 


Eq 




0.213648 


0.21473 


0.21364898 
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Table 5. Polytrope n — 1 



Quantity 


3rd order theory 


jn expansion 


CMS (N = 256) 


CMS (N = 512) 


Q 


0.089195487 


0.089195487 


0.089195487 


0.089195487 


J2 X 10^ 


1.3994099 


1.3988511 


1.3991574 


1.3989253 


-J4 X 10^ 


5.3871087 


5.3182810 


5.3203374 


5.3187997 


Je X 10^ 


3.9972442 


3.0118323 


3.0133819 


3.0122356 


-Jg X 10^ 




2.1321157 


2.1334136 


2.1324628 


Jio X 10^ 




1.7406710 


1.7418428 


1.7409925 


-J12 X 10^ 




1.5682179 


1.5693324 


1.5685327 


Jl4 X 10^ 




1.5180877 


1.5191923 


1.5184156 



